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1.  INTRODUCTION 


An  important  problem  in  turbomachinery  is  the  prediction  of 
the  flutter  boundaries  of  the  compressor  and,  in  order  to  compute 
these  boundaries,  the  unsteady  aerodynamic  forces  need  to  be 
understood.  The  aerodynamic  causes  of  flutter  in  a  compressor 
can  be  very  complex,  for  example,  the  interaction  between  the 
flows  induced  by  the  rotor  and  stator,  but  an  important  class  of 
flutter  is  caused  by  a  relatively  simple  aerodynamic  flow  due  to 
vibration  of  the  blades  on  a  given  compressor  row  without  the 
added  complexity  of  exterior  interactions.  This  latter  problem 
is  the  simplest  unsteady  aerodynamic  phenomenon  of  interest.  In 
a  compressor  the  oncoming  velocity  varies  along  the  compressor 
blade,  increasing  from  subsonic  speeds  at  the  hub  to  perhaps 
transonic  speeds  at  the  tip;  it  is  the  prediction  of  unsteady 
transonic  flow  in  a  compressor  that  is  the  most  difficult.  This 
proposal  is  concerned  with  the  prediction  of  flutter  in  a  vibrat¬ 
ing  compressor  row  when  the  flow  is  transonic. 

Because  of  the  general  complexity  of  even  the  simple  com¬ 
pressor  row  flow,  several  simplifications  must  be  made  in  order 
to  obtain  a  tractable  problem.  The  first  simplification  is  to 
assume  that  radial  velocities  are  small  relative  to  axial  veloci¬ 
ties  and  unwrap  the  compressor  row  to  form  a  two-dimensional 
cascade.  The  second  simplification  is  to  assume,  in  the  first 
instance,  that  viscous  effects  are  negligible  and  that  any  shock 
waves  present  are  sufficiently  weak  such  that  potential  theory 
can  be  used.  In  certain  circumstances  the  potential  equation 
model  can  be  simplified  further  to  yield  the  unsteady,  high  fre¬ 
quency  transonic  small  disturbance  (TSD)  equation.  This  is  the 
simplest  equation  that  will  represent  the  dominant  feature  of 
the  flow.  It  can  be  noted  that  viscous  effects  can  be  repre¬ 
sented  by  the  inclusion  of  a  boundary  layer  model. 
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A  major  problem  with  potential  theory  is  that  it  is  not 
valid  (Ref.  1)  when  shock  waves  are  present  and,  as  a  conse¬ 
quence,  that  shock  locations  are  not  predicted  adequately  for 
medium  to  strong  shock  strengths,  although  the  agreement  for 
flows  with  weak  shocks  is  acceptable.  This  disadvantage  of  the 
potential  formulation  can  be  overcome  by  modifying  the  theory 
to  give  the  correct  shock  jump.  This  has  been  done  by  Nixon 
(Ref.  2)  for  the  steady  small  disturbance  equation  and  by 
Kerlick,  Nixon,  and  Ballhaus  (Ref.  3)  for  the  unsteady  low  fre¬ 
quency  equation.  A  more  complete  analysis  for  steady  flow  (Ref. 
4)  indicates  that  such  modified  equations  can  essentially  dupli¬ 
cate  the  Euler  equations. 

Most  of  the  numerical  methods  for  predicting  unsteady  tran¬ 
sonic  flow  stem  from  the  work  of  Ballhaus  and  Goorjian  (Ref.  5), 
who  developed  an  implicit  algorithm  to  solve  the  low  frequency 
TSD  equation  for  isolated  airfoils.  An  extension  of  this  algo¬ 
rithm  to  high  frequency  motion  was  made  by  Rizzetta  and  Chin 
(Ref.  6)  .  For  cascade  flows  a  version  of  the  Rizzetta-Chin  algorithm 
was  used  by  Kerlick  and  Nixon  (Ref.  7)  to  develop  a  method  for  cas¬ 
cade  flows  for  unstaggered  cascades.  An  alternative  to  the  use 
of  the  unsteady  TSD  equation  has  been  considered  by  Verdon  and 
Caspar  (Ref.  8)  who  use  a  time  linearized  perturbation  of  the 
full  potential  equation  for  staggered  and  unstaggered  cascades 
in  a  subsonic  flow.  The  methods  of  References  5-7  can  treat 
moving  shock  waves  but  it  is  not  clear  at  present  whether  the 
time  linearized  method  is  capable  of  such  a  representation, 
since  it  does  not  possess  the  essential  double-valued  character 
at  a  shock  location. 

In  a  cascade  flutter  analysis  the  dominant  parameter  is  the 
interblade  phase  angle  which  can  determine  the  conditions  at 
which  the  blade  will  flutter.  In  a  numerical  method,  such  as 
those  discussed  earlier,  a  flutter  calculation  can  involve  com¬ 
puting  a  test  case  for  each  value  of  the  interblade  phase  angle. 
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This  is  computationally  expensive  and,  in  order  to  reduce  the 
cost,  a  technique  developed  by  Nixon  (Ref.  9)  may  be  used.  This 
technique,  which  is  based  on  indicial  theory,  decouples  the  cas¬ 
cade  problem  into  a  series  of  elementary  problems,  which  are 
parametrically  related,  and  only  one  solution  need  be  computed 
numerically.  This  elementary  solution  consists  of  a  number  of 
blades,  all  stationary  with  the  exception  of  one  blade  which 
undergoes  a  step  change.  Once  computed,  the  full  unsteady  solu¬ 
tion  for  any  frequency  and  interblade  phase  angle  can  be  con¬ 
structed  by  a  judicious  superposition. 

The  present  work  is  concerned  with  the  development  of  the 
unsteady  transonic  flow  through  a  cascade.  Earlier  work  concen¬ 
trated  on  unstaggered  cascades  with  the  present  study  concerned 
with  extending  the  method  to  treat  staggered  cascades.  These 
studies  are  performed  under  contracts  N00018-81-C-0169  and 
N00019-82-C-0288 .  In  addition  to  the  development  of  the  flow 
prediction  code,  the  ideas  of  Reference  9  have  been  tested  to 
determine  whether  flutter  prediction  at  transonic  speeds  is 
possible  without  resorting  to  an  expensive  series  of  calcula¬ 
tions  for  a  range  of  interblade  phase  angles. 

The  flow  prediction  code  is  a  development  of  the  Ballhaus- 
Goorjian  (Ref.  5)  code  to  include  the  high  frequency  algorithm 
of  Rizzetta  and  Chin  (Ref.  6)  and  to  extend  the  geometry  of  the 
problem  from  an  isolated  airfoil  to  that  of  a  cascade.  The 
present  code  is  limited  by  computer  storage  to  a  12-blade  cas¬ 
cade.  The  code  is  composed  of  two  distinct  parts,  namely  a 
steady  flow  part  which  is  periodic  in  a  blade-to-blade  sense, 
which  is  solved  by  a  relaxation  algorithm,  and  an  unsteady  flow 
part  which  uses  an  alternating  direction  implicit  (ADI)  algo¬ 
rithm.  The  unsteady  part  of  the  code  requires  that  the  inter¬ 
blade  phase  angle  be  specified  or,  if  run  in  the  indicial  mode, 
that  a  relatively  large  number  of  blades  (usually  seven)  be  used 
for  the  calculation. 
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At  present,  the  computer  code  works  for  unstaggered  cas¬ 
cades  but  the  extension  to  staggered  cascades  has  not  been  com¬ 
pleted. 

In  order  to  test  the  validity  of  the  method  of  Reference  9, 
the  flow  prediction  code  was  run  in  its  indicial  mode  and  the 
resulting  indicial  responses  used  to  compute  the  magnitude  and 
phase  lag  of  lift  and  pitching  moment  coefficients.  A  compari¬ 
son  of  values  computed  by  this  method  and  the  direct  calculation 
method  discussed  earlier  shows  good  agreement. 

The  computation  time  required  to  obtain  the  indicial 
response  for  a  seven-blade  cascade,  and  to  use  this  data  to 
compute  the  oscillatory  lift  and  moment  coefficients,  is  approxi¬ 
mately  500  cpu  seconds  on  a  CDC  7600  computer.  This  figure  is 
essentially  unchanged  if  a  range  of  interblade  phase  angles  and 
frequencies  are  to  be  computed.  The  direct  computation  requires 
250  cpu  seconds  for  each  interblade  phase  angle  and  for  each 
frequency.  Hence,  in  a  flutter  calculation  the  method  of  Refer¬ 
ence  9  offers  a  considerable  reduction  in  computation  time. 

While  it  can  be  stated  that  the  agreement  of  the  present 
methodology  with  other  available  test  cases  is  satisfactory, 
there  are  areas  that  need  to  be  further  investigated  before  the 
true  worth  of  the  technique  can  be  determined.  These  points  are 
noted  in  the  text. 

This  report  consists  of  two  more  or  less  separate  parts, 
namely  the  development  of  an  algorithm  for  the  unsteady  cascade 
code  and  the  development  of  the  indicial  theory.  Although  the 
latter  uses  the  numerical  data  generated  by  the  unsteady  cascade 
code,  it  is  not  tied  to  this  particular  program;  any  suitable 
computer  code  can  be  used  to  calculate  the  necessary  input. 
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2. 


UNSTEADY  TRANSONIC  CASCADE  CODE 


2.1  Unstaggered  Algorithm 


Although  a  preliminary  version  of  the  cascade  algorithm 
has  been  described  in  Reference  7,  there  have  been  sufficient 
modifications  to  the  algorithm  to  warrant  a  further  description. 

The  basic  equation  is  the  high  frequency  transonic  small 
disturbance  (TSD)  equation  (see  Ref.  5)  : 

Ad),  ,  +  2Bd>  ,  =  C&  +  d)  (1) 

Ytt  yxt  ^xx  Yyy 


2  2  -2/3 

where:  A  =  v  M  6  ' 

00 

B  =  Av_1 

C  =  (1  -  M2)6“2,/3  -  (y  +  l)Mq4) 

Here  <}>  is  the  perturbation  velocity  potential,  is  the  free 
stream  Mach  number,  y  is  the  ratio  of  specific  heats,  6  is  the 
thickness  parameter  and  q  is  the  transonic  scaling  parameter. 
The  quantities  x,  y,  t,  <}>  are  in  the  scaled  units  of  Reference 
5.  The  reduced  frequency  based  on  chord  c  is  v  =  uc/U^,  where 
w  is  the  frequency  and  U^  is  the  free  stream  velocity. 

Ballhaus  and  Goorjian  (Ref.  5)  solve  Equation  (1)  in  the 
case  where  A  =  0,  which  applies  for  reduced  frequencies  less 
than  about  0.2  in  the  code  LTRAN2.  The  unsteady  effects  which 
are  of  interest,  primarily  cascade  flutter,  occur  at  reduced 
frequencies  (based  on  blade  chord)  of  about  one  to  ten.  Since 
the  Ballhaus-Goor j ian  algorithm  is  not  designed  for  these  fre¬ 
quencies,  high  frequency  modifications  must  be  made. 

Houwink  and  van  der  Vooren  (Ref.  10)  and  Hessenius  and 
Goorjian  (Ref.  11)  solve  Equation  (1),  still  with  A  =  0,  but 
with  improvements  in  the  definition  of  C^,  and  the  airfoil  and 
wake  boundary  conditions.  Rizzetta  and  Chin  (Ref.  7)  solve  the 
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full  Equation  (1)  ,  including  the  A<j>  term,  with  improved 
boundary  conditions.  Their  algorithm  is  summarized  briefly 
in  the  following  section. 

2.2  Rizzetta-Chin  Algorithm 


The  Ballhaus-Goor j ian  algorithm  has  two  components,  a 
steady  TSD  equation  solver  by  either  successive  line  over¬ 
relaxation  (SLOR)  or  approximate  factorization  (AF2)  algo¬ 
rithms,  and  an  unsteady  solver  based  on  the  alternating  direc¬ 
tion  implicit  (ADI)  algorithm.  Only  the  unsteady  solver  is 
modified  for  high  reduced  frequencies.  The  airfoil  boundary 
condition  is  applied  on  the  split  coordinate  line  at  y  =  0 
(see  Fig.  1) .  The  algorithm  approximates  Equation  (1)  by  finite 

differences  at  a  time  level  tn+1//2  in  order  to  solve  for  the 

nr4“  1 

velocity  potential  function  <J>  at  time  level  t  ,  given  <f>  at 
time  tn  (Ref.  12).  For  the  high  frequency  equation,  Rizzetta 
and  Chin  approximate  the  additional  term  A<j>tt  by  a  first-order 
accurate  three-point  backward  difference: 
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n+1 


-  2 
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tt 


(At) 


+  1  At 

2 -  +  V  *ttt 


(At2) 


(2) 


n+1/2 


The  algorithm  is  divided  into  an  x-sweep  and  a  y-sweep  as  fol¬ 
lows  : 


x-sweep 


2B (At)  16x( 


<f>  —  <}>n)  .  .  =  D  f .  .  +  <{) n 

v  *  1,3  X  1,3  yy  1,3 


(3) 


y-sweep 


A  ( At )  2Un+1  -  2  <}>n  +  <)>n  1).  ■  +  2B  ( A  t)  16  (<{)n+1 

1  /  J  X 

=  1  6  Un+1  -  4,n)  .  . 

2  yy'v  i,3 
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Here,  the  term  D  f.  .is  the  nonlinear  term  defined  in 

X  1  f 

Equation  (10)  of  Reference  5  which  includes  type-dependent  dif¬ 
ferencing  by  means  of  a  switching  operator.  The  only  change 
from  the  Ballhaus-Goor j ian  method  is  the  addition  of  the  A^^ 
term  in  the  y-sweep  which  is  treated  in  such  a  way  to  retain 
conservation  form.  In  the  low-frequency  case,  this  algorithm 
ensures  correct  capture  of  a  moving  shock  wave  (Ref.  12).  The 
total  algorithm  can  be  readily  obtained  by  eliminating  the 
intermediate  variable  <f> ,  whereby  in  the  subsonic  (central  dif¬ 
ferences  6  ) 

xx 

A  (At) _2  ( <|>n+1  -  2<j,n  +  <frn_1).  •  +  2B(At)'_1((i)n+1  -  <j>n)  i ,  j 

1  r  J 


C 
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6 

xx 
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n+1 
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yy 
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n+l 


+  d)n)  .  .  +  TE 

1  r  3 


(5) 


Here,  TE  represents  the  truncation  error  incurred  by  approxi¬ 
mately  factorizing  the  differential  operators.  For  subsonic 
flow,  the  truncation  error  is  given  by 


TE 


ACAt 

4B 


6  5  ,  ,  <f> 

xx  tt 


2 

CAt 

8B 


6  16  6  6.  4> 
x  xx  yy  t 


(6) 


The  symbol  6  ^  indicates  integration.  Equation  (6)  indicates 
that  the  new  algorithm  is  first-order  time  accurate,  whereas 
the  low  frequency  algorithm  is  second-order  time  accurate. 

The  pressure  coefficient,  C  ,  is  given  by 

P 

Cp  =  -2S2/3Ux  +  v<f>t)  (7) 


and  the  airfoil  tangency  condition  is 

4>  (x,0)  =  fx  +  vft  (8) 

where  f(x,t)  represents  the  airfoil  surface.  The  downstream 
boundary  condition  is 

^x  +  v<^t  ~  0  (9) 
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and  the  initial  condition 


<j>t  (0)  =  h  (x , y )  (10) 

Usually  h(x,y)  is  put  equal  to  zero. 

The  far  field  boundary  condition  of  Ballhaus  and  Goorjian, 
namely  <j>  =  0  (reflective)  for  the  upstream  boundary  is 
applied  unchanged.  This  is  not  really  the  physically  correct 
boundary  condition  for  a  cascade,  but  the  mesh  is  stretched  to 
some  800  chords  away  from  the  airfoil,  which  should  minimize 
the  effect  of  the  boundary  condition.  This  minimum  influence 
of  the  far  field  boundaries  has  been  observed  in  numerical  exper¬ 
iments  . 

The  vorticity  shed  at  the  trailing  edge  is  represented  by  a 
discontinuity  r(x,t)  =  <J>(x,+0,t)  -  <j>(x,-0,t)  in  the  perturbation 

potential  across  the  blade  chord  line.  Following  References  6 
and  11,  this  discontinuity  is  convected  downstream  with  the  free 
stream  velocity,  thus 

r  +  vr,  =  o  (11) 

x  t 

This  ensures  that  the  pressure  is  continuous  across  the  wake. 

The  initial  steady  flow  through  the  cascade  is  calculated 
in  the  manner  described  in  Reference  7. 

The  unsteady  solution  is  periodic  over  the  whole  cascade  of 
N  blades,  so  the  unsteady  calculation  takes  place  on  a  mesh  which 
is  an  N-f old  replica  of  the  steady-state  mesh,  thus  the  entire 
cascade  is  represented.  Again,  the  blades  and  wakes  are  repre¬ 
sented  by  split  coordinate  lines. 

First,  the  grid  is  swept  in  the  y-direction  from  bottom  to 
top  with  implicit  quadridiagonal  solvers  used  on  the  y-constant 
lines.  The  lower  surfaces,  interblade  lines,  and  the  upper  sur¬ 
faces,  are  solved  in  turn. 
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Next,  the  grid  is  swept  in  the  x-direction.  For  the  region 
ahead  of  the  leading  edge,  periodic  boundary  conditions  are  used 
with  a  period  of  the  entire  grid  (N  interblade  spaces) ,  and  the 
periodic  tridiagonal  solver  is  again  used.  From  the  leading 
edge  to  trailing  edge,  each  x-constant  line  for  each  interblade 
region  is  solved  using  a  separate  tridiagonal  inversion  with 
thin  airfoil  boundary  conditions  applied  at  both  ends.  Aft  of 
the  trailing  edge,  the  periodic  solver  for  the  entire  cascade 
(N  blades)  is  used,  with  the  upper  surface  potentials  being 
obtained  by  adding  the  corrected  circulation  propagated  back 
from  the  trailing  edge  according  to  the  wake  condition,  Equation 
(11)  . 

2.3  Results  for  Unstaggered  Mesh 

A  series  of  computations  were  made  for  several  test  cases 
at  both  subcritical  and  supercritical  Mach  numbers.  In  Figure 
2  is  shown  the  steady  flow  through  an  advanced  supercritical 
cascade  designed  by  Sanz  (Ref.  13)  at  a  free  stream  Mach  number 
(M^)  of  0.692  and  a  gap/chord  ratio  (G)  of  1.373.  The  agreement 
with  the  results  of  Sanz  is  adequate.  It  should  be  noted  that 
it  proved  to  be  relatively  difficult  to  converge  this  case. 

In  Figure  3  the  unsteady  pressure  difference  over  a  flat 
plate  oscillating  in  bending  with  reduced  frequency  (v)  of  one 
and  M^  =  0.5,  G  =  1.0,  is  compared  to  the  results  of  Verdon  and 
Caspar  (Ref.  8) .  The  agreement  is  good.  For  a  similar  case 
involving  torsion,  the  agreement  between  the  two  methods  is  not 
as  good,  as  can  be  seen  in  Figure  4.  The  reason  for  this  dis¬ 
crepancy  is  not  clear  at  present.  Figure  5  shows  the  pressure 

difference  for  a  NACA0012  section,  M  =  0.5,  G  =  1.0,  v  =  1.0 

00 

oscillating  in  torsion.  Again,  the  agreement  with  the  result 
of  Verdon  and  Caspar  is  not  good. 

Figure  6  shows  the  unsteady  pressure  difference  for  a 
NACA0012  cascade  at  a  supercritical  Mach  number  of  0.75,  v  = 

1.0,  G  =  2.  Figure  7  shows  the  unsteady  pressure  difference 
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for  the  Sanz  cascade  at  M  =  0.707,  G  =  1.373,  v  =  1.0.  It  can 

oo 

be  seen  that  the  presence  of  shock  waves  affects  the  unsteady 
flow  considerably. 

A  comparison  of  the  magnitude  and  phase  for  CL  and  for 
several  cascades  is  shown  in  Table  I. 


3.  UNSTEADY  CASCADE  ALGORITHM 


In  a  staggered  cascade  the  natural  coordinate  system  is 
one  in  which  periodicity  can  be  applied  along  grid  lines.  Con¬ 
sequently,  if  coordinates  (£,n)  are  defined  as 


-1/3 

x=£+ysin0  6  ' 

■z  s 

y  =  n  cos  6g 


(12) 


where  0g  is  the  stagger  angle,  the  TSD  equation,  Equation  (1), 
can  be  written  as 


where : 


A4> 


tt 


+  2B<j> 


=  (C  + 


s2) 


-  2S 


Sn 


+  4> 


nn 


2 

A  =  A  cos  e 

s 

2 

B  =  B  cos  0 

s 

C  =  6  sin  6 

s 


(13) 


The  tangencv  boundary  condition  is  now 

(j>  U,t)  =  (fr+  vf.)  cos  0  +  S<f>  (14) 

n  t. 

Equation  (13)  is  similar  to  the  unstaggered  equation.  Equation 
(1)  ,  with  the  exception  of  the  cross-derivative,  -2S<j>^.  This 
term  can  cause  numerical  difficulties  because  of  the  fact  that 
in  the  (£,n)  coordinate  system,  the  domain  of  dependence  in  the 
supersonic  zone,  defined  by  the  characteristics,  is  no  longer 
contained  in  a  simple  difference  stencil,  as  shown  in  Figure  8. 
A  solution  to  this  difficulty  has  been  proposed  by  Rae  (Ref.  14) 
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for  steady  flow,  who  used  a  seven-point  formula  to  represent  the 
cross-derivative.  Thus, 


4 


( A  £  A  n ) 


^  di-l,  j+l+ *i  ,  j+1  +  ^i-lfj 
+  h,j-i +  t+i,j  •  ♦i+i,j-i> 


(15) 


In  the  relaxation  scheme,  data  at  the  intermediate  time  level, 
denoted  by  $ ,  is  given  by 


(ATl)’2(^i,j  +  l  "  2*ij  +  *i,j-l}  +  C(A5)  2[^i-l,j  "  2^ij  + 


-1 


■  ‘(-+i- i,j+i+t,j+i  +  h-i,j  "  2hj 


+  h,j-i  +  ♦i+i,j  '  +i+i,j-i! 


-  s2 (as)  2<i1.lij  "  2*ij  + 


(16) 


where  the  superscript  n  denotes  a  value  at  the  previous  time 
level.  The  final  result  at  time  n+1  is  given  by 


*n+1 

13 


+  ( 1-w )  ( 4>  ~  (pn) 


=  03  C 


1D 


(17) 


where  u  is  a  relaxation  factor.  In  regions  of  supersonic  flow, 
the  term  in  square  brackets  in  Equation  (16)  is  upwind  differ¬ 
enced.  None  of  the  other  £  differences  are  treated  this  way 
since  they  are  due  to  the  <j>  term  in  the  original  equation, 
Equation  (1) . 

Some  results  of  using  this  algorithm  are  shown  in  Figure  9, 
where  the  steady  pressure  distribution  through  a  cascade  of 
NACA0012  section  blades  at  M  =0.75  and  G  =  2  are  shown  for 

oo 

various  stagger  angles.  It  may  be  seen  that  the  effects  of 
stagger  are  significant. 

For  unsteady  flow  problems,  a  difficulty  arises  in  incorpo¬ 
rating  the  seven-point  scheme  into  the  alternating  direction 
algorithm  (ADI)  used  in  the  existing  code.  The  difficulty  is 
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that  the  ADI  algorithm  sweeps  first  in  the  £  direction  and  then 
in  the  n  direction,  therefore  at  a  given  time  level  all  of  the 
unknowns  in  a  row  must  be  at  the  same  time.  Thus  the  differenc¬ 
ing  of  the  cross-derivative  operator  6  ,  where 


6£n^  =  [(*i-l,j  +  *i,j+l  "  *ij  “  *i-l,j+l)/(A5An)] 


+  ^i+l,j-l  +  ^i,j-l 


>ij)/(A?An)l  (18) 


must  be  decomposed  into  the  appropriate  time  levels.  If  the 

first  term  in  square  brackets  on  the  right-hand-side  of  Equation 

SE 

(18)  is  denoted  by  <5  d>  and  the  second  term  by  6  rJ"‘<j>,  then  the  ADI 

£n  £n 

algorithm  is  as  follows: 


£-sweep 


ft  V*“*n> 


I  Scca  +  *n)  +  <Snn4>n 


n-sweep 


+  V  (i  +  *n)  -  2SAS®V 

-  S  (1-A)6  fE  ($  +  4>n) 


A  ,  n+1  0,n  , n-1.  2B  ,,n+l  j, 

— 2  ( <f>  -  2<J>  +  4>  )  +  6,(<P  “  <f>) 


At 


1  _  ,  n+1  n.  NW,  n+1  n. 

=  -6  ( 4>  -  <}>  )  -  S  A  o  ( -  9  ) 

2  nt 


(19) 


(20) 


X  is  a  parameter,  which  if  equal  to  1/2,  gives  a  second-order 
accurate  algorithm. 

On  the  upper  surface  of  the  blade  the  cross-derivative  term 
is  represented  by  6^  and  on  the  lower  surface  by  6^.  This 
leads  to  a  first-order  accurate  treatment  of  the  boundary  condi¬ 
tions  . 

At  the  present  time  the  above  algorithm  does  not  work;  there 
seems  to  be  a  problem  with  the  boundary  conditions  at  the  leading 
and  trailing  edges  of  the  blade. 
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4.  APPLICATIONS  OF  THE  INDICIAL  THEORY 


The  indicial  theory  is  described  in  Reference  9,  and  is 

given  as  the  Appendix.  In  brief,  if  an  elementary  solution 

for,  say,  the  lift  coefficient  CT  (t)  due  to  a  motion  by  one 

-Li 

blade  in  a  cascade  of  2N+1  blades°with  the  other  blades  sta¬ 
tionary,  the  lift  coefficient  on  the  blade  for  the  whole  cas¬ 
cade  moving  with  a  phase  lag  of  a  is 

N 

C_  (t)  =  l  C_  (t-io)  (21) 

L  i=-N  Li 

where  CT  (t)  is  the  lift  on  the  ith  blade  in  the  elementary 
solution.  If  the  variable  blade  undergoes  a  step  change  in  the 
elementary  solution  of  e,  then  Equation  (21)  gives  an  indicial 
response  C  (t) ,  where 

Li 

e 

CL  (t)  =  (CL(t)  “  CL  )/e  (22) 

e  o 

and  any  time-dependent  motion  can  be  constructed  using  the 
indicial  theory  of  Reference  15.  Thus,  for  an  oscillatory 
cascade  the  lift  coefficient  is  given  by 


CT  (t)  =  CT  +  CT  (t)e(o)  +  C  (t)  iLifell  dT  (23) 
L  Lo  Le  Jo  Le  dt 

where  e (t)  is  a  specified  motion.  If  the  motion  is  harmonic 
with  reduced  frequency  v,  the  real  and  imaginary  parts  of  C 
are  given  by 


ACl  (t) 

E 


sin  vxdx 


(24a) 


v 


oo 

o 


AC  (t) 


cos  vxdx 


(24b) 
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where 


ACl  (t)  =  CL  (»)  -  CL  (t)  (25) 

£  £  £ 

R  I 

A  similar  relation  appears  for  the  moment  coefficient,  C^,  C^. 

Turbomachinery  flutter  will  occur  in  simple  torsion  if 
is  positive  and  for  simple  bending  if  is  positive. 

J_i 

4.1  Results 

The  above  theory  was  used  to  compute  lift  and  moment  coeffi¬ 
cients  for  cascades  with  variations  of  interblade  phase  angle, 
frequency  and  Mach  number.  The  elementary  solutions  necessary 
for  the  evaluation  of  the  sum  in  Equation  (21)  are  computed 
using  the  unstaggered  cascade  code  described  earlier. 

In  order  to  test  the  method,  the  indicial  theory  was  used 
to  compute  the  magnitude  and  phase  angle  for  a  cascade  of 

NACA0012  sections  in  torsion  at  v  =  1.0,  G=2  and  at  M  =0.5 

00 

and  0.75.  These  were  compared  to  directly  computed  results,  and 
this  is  shown  in  Table  2. 

The  variation  of  the  real  and  imaginary  moment  coefficients 
with  interblade  phase  angle  for  a  NACA0012  cascade  at  =  0.75 
and  a  flat  plate  cascade  at  =  0.5  are  shown  in  Figure  10.  It 
can  be  seen  that  there  is  a  considerable  difference  between  the 
results . 

Figure  11  shows  the  variation  of  real  and  imaginary  moment 
coefficients  for  a  NACA0012  cascade  at  =  0.75,  G  =  2  for  vari¬ 
ous  frequencies.  It  can  be  seen  that  the  variation  with  frequency 
is  considerable  and  that  at  v  =  1.2  the  imaginary  part  of  the 
pitching  moment  becomes  positive,  indicating  a  tendency  to  flutter. 

It  proved  difficult  to  obtain  indicial  responses  to  the 
supercritical  cascade.  The  reason  is  not  clear,  although  it  is 
possible  that  one  of  the  nonunique  solutions  described  by 
Steinhoff  and  Jameson  (Ref.  16)  may  be  present. 
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The  computation  time  required  to  obtain  the  indicial 
response  for  a  seven-blade  cascade,  and  to  use  this  data  to 
compute  the  oscillatory  lift  and  moment  coefficients,  is  approx¬ 
imately  500  cpu  seconds  on  a  CDC  7600  computer.  This  figure  is 
essentially  unchanged  if  a  range  of  interblade  phase  angles  and 
frequencies  are  to  be  computed.  The  direct  computation  requires 
250  cpu  seconds  for  each  interblade  phase  angle  and  for  each 
frequency.  Hence,  in  a  flutter  calculation  the  method  of  Refer¬ 
ence  9  offers  a  considerable  reduction  in  computation  time. 


5.  CONCLUDING  REMARKS 

A  computer  code  to  calculate  the  unsteady  transonic  flow 
through  a  cascade  has  been  developed  for  unstaggered  cascades. 
Although  an  algorithm  for  staggered  cascades  has  been  derived, 
it  has  not  yet  resulted  in  a  workable  computer  code.  Results 
for  the  unstaggered  case  agree  quite  well  with  other  solutions 
for  subsonic  flows,  but  there  are  some  discrepancies  that  should 
be  explained  in  the  future. 

This  computer  code  has  been  used  to  generate  input  data  for 
testing  the  cascade  indicial  method  of  Reference  9.  The  results 
look  encouraging,  but  more  testing  is  desirable.  The  indicial 
method  promises  greatly  reduced  computer  time  for  cascade  flutter 
calculations . 
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TABLE  1.-  CL  AND  CM  FOR  VARIOUS  CASCADES 


Section 

M 

oo 

G 

V 

CL 

^L 

CM 

M 

^M 

Flat  Plate 

0.7 

2.0 

1.0 

8.61 

-26.2 

2.068 

-63.3 

NACA  0012 

0.75 

2.0 

1.0 

3.17 

-57.5 

1.648 

-231.4 

Sanz 

0.707 

1.373 

1.0 

2.90 

-51.6 

1.132 

-221.5 

TABLE  2.-  COMPARISON  OF  MAGNITUDES  AND  PHASE  ANGLES  OF 
LIFT  AND  MOMENT  COEFFICIENTS  FOR  NACA  0012  SECTION 
OF  INTERBLADE  PHASE  ANGLE  OF  180° 


M 

OO 

Method  of  Ref.  1 

Direct  Calculation 

CL 

i  o 
^L 

Cw 

M 

rO 

^M 

CL 

O  1-q 
-e- 

CM 

-O 

^M 

0.50 

6.840 

-13.2 

2  .040 

-173.3 

6.776 

-  9.4 

1.968 

-16  9 . 7 

0.75 

3.428 

-57.7 

1.812 

-227 .6 

3.170 

-57.5 

1.647 

-231.4 
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Figure  1.  LTRAN2  mesh 
conditions  for  an 


schematic  and  boundary 
isolated  airfoil. 
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Figure  2 


Steady  pressure  distribution  around  a 
supercritical  blade  cascade; 

Moo  =  0.692  ,  G  =  1.372. 
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Unsteady  pressure  difference  coefficient,  Ap 


Flat  plate  bending 


Figure  3.-  Pressure  distribution  for  a  flat  plate 
oscillating  in  bending;  Ma  =  0.5,  v  =  1.0. 
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Unsteady  pressure  difference  coefficient 


Figure  4.-  Unsteady  pressure  distribution 
flat  plate  in  torsion;  =  0.5,  G  =  1.0, 


around  a 

v  =  1.0. 
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Figure  5.-  Unsteady  pressure  distribution  around  a 
cascade  of  NACA  0012  blade  in  torsion: 

=  0.5,  G  =  1.0,  v  =  1.0. 
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Distance  from  leading  edge,  x 

NACA  airfoil  0012,  torsion 

Pressure  distribution  for  a  NACA  0012  blade 
Hating  in  torsion;  Mro  =  0.75,  v  =  1.0. 
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Unsteady  pressure  difference  coefficient  Ap 


Figure  7.-  Unsteady  pressure  distribution  around  a 
supercritical  blade  cascade  in  torsion; 

Mo,  =  0.70  7,  G  =  1.3  72  ,  v  =  1.0. 
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k-I 


Figure  8.-  Schematic  of  difference  cell  in  a 
skewed  coordinate  system. 
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Pressure  Coefficient 
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CMR 


degrees 


Interblade  Phase 
Angle 


(a)  Real  part 

Figure  10.-  Variation  of  pitching  moment  with 
interblade  phase  angle . 
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Figure  10.-  Concluded. 
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Figure  11.-  Variations  of  moment  coefficients  with  interblade  phase  angle 
and  frequency.  NACA  0012  section  Moo  =  0.75,  G  =  2.0. 


o  O 


^uauioui  BuTtp^xd  jo  :ppd 


■J 0  O 

■  i 

i—l  IN 

I  I 

AjeutShuii 


33 


Figure  11.-  Concluded. 


APPENDIX 


NOTES  ON  THE  COMPUTATIONS  OF  UNSTEADY 
TRANSONIC  CASCADE  FLOWS 
David  Nixon* 

Nielsen  Engineering  &  Research,  Inc. 

Mountain  View,  California 

INTRODUCTION 

In  an  oscillating  cascade  there  is  by  definition 
a  fundamental  periodicity  that  occurs  every  2N  blades 
(2N  is  the  number  of  blades  in  the  compressor  row) . 
The  unsteady  flow  at  each  blade  will  have  a  periodic 
boundary  condition,  as  in  steady  flow,  but  will  lag 
by  a  phase  angle  of  (p/N)ir  in  relation  to  the  neigh¬ 
boring  blade,  where  p  is  an  integer  less  than  or 
equal  to  N  whose  value  is  determined  as  part  of  the 
flutter  solution.  In  a  nonlinear  transonic  numerical 
scheme  the  choice  is  between  computing  the  entire'  2N 
blade  sequence  with  the  usual  periodic  boundary  con¬ 
ditions  at  the  extremities  or  to  compute  the  usual 
three  blade  cascade  problem  for  each  specified  value 
of  p.  These  numerical  calculations  are  computation¬ 
ally  expensive  and  it  is  desirable  to  reduce  the 
overall  cost  of  a  flutter  calculation.  Both  of  these 
choices  involve  a  large  amount  of  computer  time  for 
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practical  cases,  and  in  case  of  the  first  choice,  a 
major  development  of  a  computer  code.  However,  it 
is  possible  to  devise  a  simpler  approach  to  the 
problem. 

The  basic  idea  of  this  paper  is  to  devise  an  ele¬ 
mentary  problem  in  which  only  one  blade  is  in  motion, 
the  others  being  fixed;  the  motion  may  be  any  general 
time  dependent  function.  This  removes  the  problem  of 
computing  the  flow  for  each  blade  phase  angle.  This 
elementary  problem  is  solved  for  a  particular  moving 
blade  and  the  functional  form  of  the  velocity  poten¬ 
tial  for  both  space  and  time  is  then  known.  Because 
of  periodicity  in  both  space  and  time,  these  elemen¬ 
tary  solutions  can  be  superimposed,  with  reparameter¬ 
ized  space  and  time  variables,  so  that  the  sum  of  the 
solutions  satisfies  both  the  basic  differential  equa¬ 
tion  and  the  correct  boundary  conditions  on  every 
blade.  Although  the  most  important  application  of 
this  superposition  principle  is  the  development  of 
the  correct  linearization  of  the  transonic  flow  equa¬ 
tions  with  discontinuities,  these  equations  are  much 
too  complicated  for  an  illustration  of  the  superposi¬ 
tion  technique.  Hence,  in  the  following  only  a  simple 
subsonic  problem  is  examined.  The  general  theory  is 
directly  extendable  to  discontinuous  transonic  flows 
using  the  strained  coordinate  theory  of  Nixon  (ref.  1) 
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ANALYSIS 


Consider  the  cascade  of  2N  blades  pictured  in 
fig.  1,  where  blade  J+N  and  blade  J-N  are  identified 
The  equation  for  the  perturbation  velocity  potential 
<p(xrt)r  due  to  a  time  dependent  behavior  of  the 
blades,  will  be  linear  and  can  be  represented  as 
follows: 


L  (x,t)  <j>  =  0 


(1) 


Here,  x  is  a  general  vector  coordinate  centered  on 
the  blade  of  interest,  t  is  time,  and  L(x,t)  is  a 
differential  operator  with  variable  coefficients. 
These  coefficients  arise  from  the  mean  steady  flow 
solution  about  which  the  flow  is  perturbed,  are 
functions  only  of  x,  and  are  periodic  in  space  with 
period  of  2N  blades.  The  time  t  only  occurs  through 
the  inclusion  of  time  derivatives.  Such  an  operator 
occurs  in  the  subsonic  flow  analysis  of  Verdon  and 
Caspar  (ref.  2)  who  use  the  equation 
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where  0  is  the  steady  state  potential.  For  a  stag¬ 
gered  cascade  the  upstream  and  downstream  boundary 
conditions  will  depend  on  the 
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particular  region  of  the  flow  that  is  being  consid¬ 
ered,  .for  example  the  far  upstream  boundary  D  ,  and 

U  J 


the  far  downstream  boundary  ^  shown  in  figure  1. 

The  boundary  conditions  will  depend  on  the  steady 
far  field  velocity,  which  has  spatial  periodicity 
and  on  the  phase  of  the  waves,  produced  by  the  oscil¬ 
lating  cascade.  Since  the  individual  blade  flows 
are  identical  except  for  phase  lag  of  jt  ,  it  follows 
that  the  far  field  boundary  conditions  are 
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where  x  ,  x_  denote  the  location  of  the  upstream 
u  D 

and  downstream  boundaries  respectively  and  g^  and 
g^  are  the  corresponding  spatially  asymptotic  values 
of  £(x,t).  The  time  t.  is  given  by  =  t  -  jtQ 

where  t  is  the  interblade  phase  lag.  The  tangential 
boundary  conditions  are 


v .  (±x,  t) 


f  .  (±x ,  t . )  on  Sj  ;  j  =  -N,N  (3) 


where  v.  is  the  normal  velocity  component,  f ^  is  a 
function  of  the  specified  perturbation  of  blade 
geometry  and  denotes  the  location  of  the  blade 

surface.  The  ±  signs  denote  conditions  on  the  upper 
and  lower  surface  of  the  blades,  respectively.  The 
effect  of  the  interblade  phase  lag  is  seen  in 
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equations  (2)  and  (3)  where  the  potentials  <j:  ,  4>D  or 
the  velocity  of  the  j  *  blade  at  time  t  are  given  in 
terms  of  a  function  at  a  time  t_.  .  For  blade  number  0 
periodicity  can  be  applied  on  j  =  ±N.  In  addition  to 
the  above  boundary  conditions  the  wake  condition  is 


r 


AC  (x,t)  =  0"bn  w. 
P  3 


(4) 


where  AC^  is  the  pressure  jump  across  the  wake  of  the 

jth  blade,  w ^ ,  and  is  a  function  of  cf>(x,t).  Now  if 

the  principle  of  superposition  holds  then  a  variable 

d>.(x.,t)  can  be  introduced  such  that 
3  3 


N 

4>  (x,  t)  =  ^  j  <xj  > 

j=-N 


(5) 


whe 


X  .  =  X  -  X  . 

3  3 


and  x.  is  the  location  of  the  general  coordinate 

i  _ 

system  centered  on  the  j  blade  and  xq  =  x.  Assume 

that  6 .  is  such  that 
3 


L  (x,  t)  4)  j  =  0 


(6) 


It  is  obvious  that  if  equation  (6)  is  used  with 
equation  (5)  then  equation  (1)  is  recovered. 

The  operator  L(x,t)  is  a  dif rerential  operator 
with  coefficients  that  are  functions  of  the  steady 
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flow  over  a  cascade  and  hence  L(x,t)  is  periodic  in 
the  spatial  dimensions.  Thus  L(x,t)  is  the  same  no 
matter  on  which  blade  the  coordinate  system  is 
centered.  Hence 

L  (x ,  t )  =  L(.Xj,t)  (7) 

Furthermore,  since  t  only  appears  in  the  operator 
through  a  derivative  then  the  variable  t  can  be  re¬ 
placed  by  tj .  Thus 

L  (x,  t)  =  L  (x_.  ,  t  j  )  (8) 

Using  equations  (6)  and  (8)  then  gives 

L(x.,t .)  <|>.  =  0  (9) 

3  .3  3 

The  next  task  is  to  determine  suitable  boundary 
conditions.  Let 

‘t'-jlVV  =  9u(xu,t3)  on  Duj 

VVV  =  V*D'V  on  % 

(10) 

<$>■  (xurtj)  =  0  on  Duk,  k  ?  j 
4) .  (xD,t_.)  =  0  on  DDk,  k  jt  j 
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vj(±x.,t.)  =  y^(±Xj,t.)  on  S. 


V  . 
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v  . 
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±<Xj  -  xk),tj 


±(Xj  +  X£),tj 


=  0  on  S  .  k  =  1 ,  N  -  1 

j+k 


=  0  on  l  =  1 / N  -  1 
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AC  .  (x  .  ,  t  . )  =  0  on  w . 
PD  3  3  D 


where 


N 

cP(S,t>  =  C  cpj(h,t)  (12) 

-N 

It  should  be  noted  that  in  the  (x^,t^)  coordinate 
system  each  of  elementary  problems  for  ^  is  identical. 

Periodicity  is  applied  in  the  j  =  ±N  blades  and 
their  wakes.  The  wake  boundary  condition  is  given 
by 

AC  | (x.  -  xk) ,tj |  =  0  on  w.+k  k  =  1,N  -  1 

(13) 

ACp  | (Xj  +  x£) ,tj i  =  0  on  Wj_£  L  -  1,N  =  1 

This  is  equivalent  to  keeping  all  blades  stationary 

except  the  jth  blade  and  allows  the  relevant  wave 

transmission  through  each  blade  wake.  Note  that  the 

t  h 

time  is  always  t . ,  the  time  associated  with  the  j 
blade.  No  interblade  phase  lag  is  required  at  this 
stage . 
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conditions,  the  problem  defined  by  equations  (5), 
(10),  (11),  (12),  and  (13)  is  identical  to  the  prob¬ 
lem  defined  by  equations  (1) ,  (2) ,  (3) ,  and  (4) . 

Hence  the  problem  for  any  time  lag  t  between  blades 
can  be  constructed  from  the  _superposition  of  the 
elementary  problem  defined  by  equations  (9) ,  (11) , 

(12) ,  and  (13) .  The  superposition  mechanism  is  as 
follows . 

Let  the  solution  to  the  elementary  problem  for 
j  =  0  be  given  by  4>Q(x,t).  The  solution  for  <j>.  (x,tj) 
is  then  given  by 

<f>j(Xj,tj)  =  tQ(x,t)  (14) 

Since  the  functional  form  of  with  both  x  and  t  is 
known  from  the  elementary  solution  for  j  =  0  this 
reparameterization  is  trivial.  The  final  solution 
for  the  zeroth  blade  is  then  given  by  equations  (5) 
and  (14);  thus 

N 

4>(x,t)  =  4>0(x  -  x.,  t  -  jtQ)  (15) 

-N 

Thus  the  complete  time  dependent  cascade  flow  for 
any  interblade  phase  angle  can  be  constructed  by 
superposition  of  the  elementary  problem  defined  by 
equations  (5),  (10),  (11),  (12),  and  (13). 
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If  4>(x,t)  and  its  derivatives  are  continuous  a 

similar  relation  to  equation  (15)  can  be  constructed 

for  the  pressure  coefficient  C  (x,t).  In  addition, 

I? 

similar  formulae  can  be  derived  for  the  lift  and 
moment  coefficients. 

The  idea  discussed  above  can  be  extended  to  dis¬ 
continuous  transonic  flows  using  the  method  of 
strained  coordinates  (ref.  1) .  Results  are  given 
by  Kerlick  and  Nixon  (ref.  3) . 

CONCLUDING  REMARKS 

The  main  contribution  to  the  computation  time  for 
an  unsteady  calculation  of  cascade  flutter  in  a 
transonic  flow  is  the  need  to  repeat  the  calculation 
for  a  range  of  interblade  phase  angles.  The  present 
analysis  shows  how  this  problem  can  be  eliminated  by 
a  judicious  choice  of  elementary  solutions. 
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Figure  1.-  Sketch  of  Cascade  Problem 
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